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The "Scalar Timing Law," which is a temporal domain generalization of the well known 
Weber Law, states that the errors estimating temporal intervals scale linearly with the 
durations of the intervals. Linear scaling has been studied extensively in human and animal 
models and holds over several orders of magnitude, though to date there is no agreed 
upon explanation for its physiological basis. Starting from the assumption that behavioral 
variability stems from neural variability, this work shows how to derive firing rate functions 
that are consistent with scalar timing. We show that firing rate functions with a log-power 
form, and a set of parameters that depend on spike count statistics, can account for 
scalar timing. Our derivation depends on a linear approximation, but we use simulations 
to validate the theory and show that log-power firing rate functions result in scalar timing 
over a large range of times and parameters. Simulation results match the predictions of 
our model, though our initial formulation results in a slight bias toward overestimation that 
can be corrected using a simple iterative approach to learn a decision threshold. 

Keywords: scalar timing, Weber's law, temporal intervals, temporal coding, neural dynamics 



1. INTRODUCTION 

Errors estimating the intensity of a stimulus commonly scale lin- 
early with the magnitude of the stimulus. This relationship, called 
Weber's Law, has proven to be a surprisingly general property 
of the brain that accurately describes perception across sensory 
modalities (Weber, 1843; Coren et al., 1984). We have previously 
used basic principles to argue that this scaling naturally emerges 
if neural processes representing stimulus magnitudes have tuning 
curves with a specific mathematical form and that the generality 
of the law implies that this is a fundamental organizing principal 
of neural computation (Shouval et al., 2013). 

An analog of Weber's law in the temporal domain, called lin- 
ear scaling or scalar timing, states that errors estimating temporal 
intervals scale linearly with the duration of the intervals (Gibbon, 
1977; Church, 2003). Temporal perception has been extensively 
studied by psychologists and neuroscientists for over 150 years, 
starting in the 1860s with Fechner (Fechner, 1966), leading to 
considerable knowledge about the behavioral aspects of temporal 
perception. Much less is known, however, about the underly- 
ing neural substrate responsible for engendering observed timing 
behavior. 

Over the years many theories have been proposed to account 
for scalar timing. The scalar expectancy theory (Gibbon, 1977; 
Church, 2003) is based on a counter and accumulator model, con- 
ceptually similar to counting the ticks of a mechanical clock, and 
variability arises from comparison errors with remembered ref- 
erence values. Another class of models assumes an ensemble of 
neurons oscillating at different frequencies, and timing is pro- 
duced by decision neurons which become active only when a 
precise set of the oscillating neurons are coactive (Matell and 



Meek, 2000). These models are akin to a Fourier transform of 
the desired temporal response profile. Variability stems from the 
addition of stochastic noise to the ensemble dynamics, and it 
has recently been shown analytically that general addition of 
noise at various levels in the model can result in scalar tim- 
ing (Oprisan and Buhusi, 2011, 2014). Note that these models 
are currently derived based on continuous dynamical systems, 
not spiking neural models. Drift diffusion models have also 
been proposed to provide a mechanistic basis of interval timing, 
though with spike-statistics that are inconsistent with scalar tim- 
ing. Recent derivations of this model, with drift that is driven 
by opponent inhibitory and excitatory processes, can account 
for scalar timing (Rivest and Bengio, 2011; Simen et al., 2011). 
A final class of models, including this work, assume that tim- 
ing is derived from the state of dynamic neural responses. For 
example, time can be estimated from the threshold crossing of 
decaying neural response (Staddon et al., 1999), or from a pre- 
cisely designed set of leaky integrators (Shankar and Howard, 
2012). 

Here, by extending our earlier analysis (Shouval et al., 2013) to 
the temporal domain, we explore the relationship between neural 
dynamics and temporal perception and propose a theory of scalar 
timing based on experimentally verifiable physiological processes. 
Our approach is based on the assumption that estimates of a tem- 
poral interval vary on a trial-by-trial basis due to spike count 
variability (Dean, 1981; Tolhurst et al, 1983; Churchland et al., 
2010). Although our analysis is mathematically homologous to 
the intensity variable case, the physiological substrates of time and 
intensity estimation are quite different. Here we show that neural 
processes with with activity levels that dynamically progress with 



Frontiers in Human Neuroscience 



www.frontiersin.org 



June 2014 | Volume 8 | Article 438 | 1 



Shouval et al. 



Scalar timing and neural dynamics 



a log-power temporal profiles can account for scalar timing, in 
much the same way that our previous work showed that Weber's 
law results from log-power tuning curves. Our analysis is based on 
a linear approximation, with results that are less precise than we 
found when analyzing intensity variables. Though the log-power 
model does produces scalar timing, our initial formulation of the 
model with the linear approximation results in a small bias toward 
underestimation and slightly less variability than found in sim- 
ulations. In this paper we also derive a discrete approximation, 
which is applicable only in the temporal domain, that precisely 
estimates simulated neural variability. We also demonstrate that 
the bias is a consequence of our initial threshold selection cri- 
teria that can be easily eliminated with a simple algorithm that 
learns the correct threshold value to accurately decode desired 
intervals. 



2. METHODS AND RESULTS 

We start with the assumption that the brain uses the tem- 
poral evolution of neural activity, which progresses with pre- 
dictable stochastic dynamics, to estimate intervals. Specifically, 
we assume that some stimulus at time t = 0 initiates a neural 
process (describing either a single neuron or, more likely, a neu- 
ral ensemble) that is characterized by a spike rate function r(t) 
(see Figures 1B,C), which either increases (Roitman and Shadlen, 
2002) or decreases (Shuler and Bear, 2006) monotonically. The 
average spike count within a window r is: 



a ■ t. 



(3) 



f t+T/2 

R(t) = / r{t+t')dt' 

Jt-x/2 



(1) 



which can be approximated as R(t) ~ rr(f) for small r values. 
As illustrated in Figure 2, the temporal interval described by 
this process is defined as the time required for R(t) to reach 
some threshold Rq, which can be set as Rq = r(t tar ) for a target 
time t tar . Noise driven fluctuations of r(t) result in variability of 
trial-by-trial estimates, t est , of the encoded target time. 

In this framework there is a direct relationship between 
the magnitude of temporal estimate errors (which can be eas- 
ily recorded using standard psychophysical methods) and spike 
count statistics that can be used to infer a mathematical form of 
r(f ), and thus the underlying physiology, subject to the linear con- 
straint on estimate errors as a function of t tar specified by the 
scalar timing law. A simple linear approximation (illustrated in 
Figure 1A) of this relationship between the interval estimate and 
spike count has the form: 



attest) 



OR{t t ar) 0 R (t tar ) 



\R'(ttar)\ \rr'(t tar )\ 



(2) 



where a#(f) is the standard deviation of the average spike count 
at time f, and R'(t tar ) is derivative of the spike count curve with 
respect to the time, estimated at the target time. Note that a t {t est ) 
is the standard deviation of the estimation of the time f over many 
trials. 

The scalar law states that errors estimating f scale linearly with 
t . Using standard deviation as the error measure: 



where a specifies the slope of linear scaling, equivalent to the 
"Weber fraction." 

Combining Equations 2 and 3: 



dr{t) dR(t) a R {t) 
r — - — = — - — = ±- 



dt 



dl 



at 



(4) 



where the + sign is valid when the slope of r(f) is positive and 
the — sign when it is negative. 

We assume that spike count variability can be characterized 
using a power-law model with the form: 



a R (t) = p (rr(t))' 



(5) 



where the parameters and p specify the specific noise model. 
This power-law model can account for many forms of spike-count 
variability. For example, p = 1/2 and yS = 1 result in Poisson 
noise, and the p = 0 case is the constant noise case, which means 
spike count variability does not depend on the spike count. 
Experimentally spike count variability is found to be close to 
Poisson and often with somewhat larger variability than Poisson 
(p ~> 1/2). Although the power-law noise is a relatively general 
model, it obviously cannot account for all forms of noise. 

Applying this form to Equation 4, we obtain a differential 
equation relating the neural firing rate to specific noise and 
estimate error models: 



P- l \ r{t)P 



dr(t) 
dt 



The solution of Equation 6 has a log-power form: 
r(t) = K ■ (±log(t/t 0 )) n 



(6) 



(7) 



where K = — p)/a) n , and n= j^-. This relationship 

holds whether r(t) rises (f > fo, "+"case) or falls (f < to,"— " 
case) monotonically. The integration constant fo has a simple 
interpretation: it is the minimal (or maximal) time that can 
be estimated using this specific monotonically increasing (or 
decreasing) firing rate function (as shown in Figures 1B,C) for 
different values of fo. Note that all the parameters of the log- 
power function are determined by measurable spike statistics and 
behavioral performance; none of them are free parameters. 

The specific shape of the general log-power form depends pri- 
marily on the spike count statistics. In the constant noise case 
(p = 0) this equation reduces to Fechner's law (Fechner, 1966). 
Hence, Fechner's law can be seen as making an implicit constant 
noise assumption. In the special and unrealistic case of propor- 
tional noise (p = 1) a power-law solution is obtained (Stevens, 
1961). 

Experimental recordings are often characterized by a nearly- 
linear relationship between mean spike count and variance(Dean, 
1981; Tolhurst et al., 1983; Churchland et al, 2010). In the near- 
Poisson case, (p = 1/2), the log-power form has an exponent 
of n = 2 (note, the examples in Figures 1B,C assume Poisson 
statistics). We have previously shown for the case of magnitude 



Frontiers in Human Neuroscience 



www.frontiersin.org 



June 2014 | Volume 8 | Article 438 | 2 



Shouval et al. 



Scalar timing and neural dynamics 




FIGURE 1 | Scalar timing and neural statistics. (A) A local linear 
approximation (green line. Equation 2) of the the average firing rate ff(f) (real 
distribution shown schematically by the gradient as a function of f, mean and 
standard deviations indicated by dashed-white and solid red lines) together 
with the scalar timing law leads to Equation 4, the solution of which (Equation 



7 for the case of Poisson noise) is the firing rate curve r(f). Note, Pf is the 
slope of the linear approximation to R(f). (B,C) Example firing rate curves 
with Poisson spike statistics for different values of the integration constant 
fo. (B) Increasing solutions are defined above minimal values at to. (C) 
Decreasing solutions are defined below maximal values at fo. 




FIGURE 2 | Temporal interval estimation. (A) A stimulus (at time t = 0) 
initiates a neural process with a mean firing rate (black line, determined by 
linear approximation theory) that decreases with time. In each trial the actual 
number of spikes varies stochastically; three trial-by-trial examples of the 
spike count variable are shown by the colored lines. The time estimate in 
each trial is determined by the first threshold crossing (Ry - horizontal dashed 
line) of the spike count variable. The actual estimated time for one trial (f est ) 
is shown in comparison to the target time (t tar ). (B) The mean time predicted 



by the model ({t as t), averaged over 200 trials) as a function of the target time. 
Blue circles based on simulations, red circles using discrete approximation. 
(C) The standard deviation of the time estimate (<jj) as a function of the 
mean predicted interval. (D) When rescaled by the mean estimated time 
(values specified by the color-code shown in the legend), the cumulative 
distributions of the actual response times overlap and are statistically 
indistinguishable (KS-test). These distributions were generated using Poisson 
statistics, a decreasing log-power function, to = 10 and r = 0.1 sec. 
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estimation that Weber's law can be based on the tuning curves of 
either single neurons or neural ensembles (Shouval et al, 2013). 
Similarly here, while scalar timing can result if a single neuron's 
time-varying activity follows a log-power function, it is more likely 
to arise from the combined activity of a heterogenous population 
of neurons whose collective activity has the appropriate form. 

To test the validity of our theory, we simulated a stochastic 
neural process with a monotonically falling spike rate in time 
(the increasing case, not discussed, is similar). Specifically, as 
per the derivations above, simulations were performed by gen- 
erating spikes using a non-homogeneous Poisson process with 
a firing-rate parameter that decreased as a log-power function 
of time. Firing rates were determined by convolving the resul- 
tant spikes trains with a square window of width r = 100 ms. 
The estimate of the temporal interval (t est ) was defined, on a 
per-trial basis, as the time at which the firing rate first reached 
threshold (-Rrj set to r{t tar )). The result of these simulations are 
shown in Figures 2B,C. The mean value of t est is close to, but 
a bit shorter than that predicted by the linear approximation 
theory (Figure 2B, blue circle) and the standard deviation is 
a linear function of the mean estimated time (Figure 2C-blue 
circle), although with a slope lower than that predicted by 
the linear approximation. Nevertheless, the rescaled distribu- 
tions are almost completely overlapping (Figure 2D) and paired 
Kolmogorov-Smirnov tests find that the differences between these 
distributions are not statistically significant (although small dif- 
ferences may emerge with longer time spans, larger a values, 
or more trials). This shows that the log-power firing rate func- 
tion indeed produces scalar timing, but that the theory described 
above results in an overestimate of error and a small bias of 
the mean. 

It is possible to obtain a discrete approximation that better cap- 
tures the simulation results. This approximation is obtained by 
dividing time into non-overlapping time bins of length r, such the 
average spike count within a time bin (designated by the integer i) 
is -R; ^ r • r((i— 0.5) • t). Under the assumptions that the bins 
are non-overlapping and have no significant temporal correla- 
tions, the spike counts in each bin are conditionally independent. 
Then, for a given spike generation model, the probability of a 
threshold crossing within a time bin as time unfolds is: 

Pc(i\Ro) = Yl p s( n \ R i)> (») 

n<R 0 

where P s (n\R) is the probability of emitting n spikes given the 
mean spike count R (note that this formulation assumes a 
decreasing function r(t), for the increasing function the sum is 
over n > Rq). The probability that the first threshold crossing 
occurs in time bin is 

P FC (j\Ro) = PSRoMzld - Pcd\Ro))- (9) 

This distribution can be used to calculate the mean, (t est ), and 
standard deviation, aj, of elapsed time estimates. Results of these 
calculations, (Figures 2B,C, red circles) agree closely with the 
numerical simulation results. The small discrepancy between this 
calculation and the actual simulations arises from partitioning 



time into non-overlapping time bins. Coarse grained simulations 
in which zero crossings are allowed only at these discrete points 
agree perfectly with the results of this discrete approximation. 

Despite the close agreement between theory and simulations, 
the model as described consistently underestimates the mean tar- 
get time. This bias, which results from the somewhat arbitrary 
decision to select Rq = r{t tar ), can be corrected if the thresh- 
olds are learned rather than chosen directly from the spike count 
function. To learn the threshold (Rq) we used a simple iterative 
learning rule: 

= ±r,{t tar = t^) (10) 

at 

where the + sign corresponds to the monotonically falling firing 
rate case, and the — sign is used for the monotonically increas- 
ing cases, and r\ < < 1 is the learning rate. This procedure quickly 
converges to provide an unbiased estimate of the target times 
(Figure 3A); error still scales linearly with time (Figure 3B) and 
the discrete approximation accounts well for the slope. 

3. DISCUSSION 

Relating behavior to its underlying physiological mechanism is a 
fundamental aim of neuroscience. Here we have shown how, in 
a class of models based on the idea that time is estimated based 
on the dynamic state of neural processes, to relate scalar timing 
to the time varying firing rate of neurons. We show that, given 
firing rate statistics characterized by a power-law, scalar timing 
arises from a log-power firing rate function with parameters that 
depend on the spike statistics. Our derivation relies on a linear 
approximation, but we have also shown that a log-power function 
results in scalar timing irrespective of this approximation. The 
initial method for setting the detection threshold using the mean 
of the firing rate function causes a small estimate bias, but this 
can be corrected using an iterative procedure to find an appro- 
priate detection threshold. Further, we have shown how to use a 
discrete approximation to calculate better estimates of encoded 
time and variability given a log-power firing rate function. These 
result depend on a mathematical analysis which is similar to the 
one used in our previous analysis of Weber's law for intensity vari- 
ables (Shouval et al, 2013). Though it produces scalar timing, 
the linear approximation is less precise in the temporal domain 
than it was when we used it to analyze the coding of stimulus 
intensity. Accordingly, it was necessary to introduced a discrete 
approximation in order to accurately calculate simulated variabil- 
ity. We also showed how the decision-selection threshold can bias 
the encoded interval. Most models that account for Weber's law 
in the intensity domain are completely distinct from models that 
account for scalar timing. Our results show that a single mathe- 
matical approach can provide a unified explanation for these two 
distinct observations. 

Although we derived firing-rate functions for the case of 
perfectly-linear scalar timing, the same procedure used to gener- 
ate Equation 4 can also be used when if the relationship between 
estimation error and time is non-linear (Grondin, 2012). Similar 
derivations are possible for other functional forms of Equation 3. 
Indeed there are experiments showing that scalar-timing is pre- 
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cisely linear only in a limited range (Getty, 1975; Bizo et al., 
2006; Grondin, 2012), and the exact forms of scaling observed 
in these experiments could be used to replace the linear scaling 
assumed here. There is no guarantee that an analytical result can 
be derived in such cases, but numerical solutions are always pos- 
sible. Similarly, the same type of approach can be used to obtain 
a firing rate function, either analytically or numerically, assuming 
non-power-law forms for neural spike statistics. 

Our analytical derivation produces a monotonically falling 
functions that are valid only below and upper threshold (to) 
and monotonically increasing functions valid only above a lower 
threshold (to). There is a lower threshold below which we can 
not evaluate time intervals, which would indicate that the mono- 
tonically increasing results are possibly more realistic. However, 
experimental results showing different Weber fractions at differ- 
ent temporal intervals could also be interpreted as indications that 
different processes are used for different time scales (Getty, 1975). 
One possibility is that falling functions could be used for very 
short temporal durations, on the order of a second or less, and 
increasing functions for longer durations. 

As we outlined above, various models of interval timing 
have been proposed over the years to account for scaler timing 
(Gibbon, 1977; Matell and Meek, 2000; Church, 2003; Durstewitz, 
2003; Oprisan and Buhusi, 2011, 2014; Rivest and Bengio, 2011; 
Simen et al., 2011; Shankar and Howard, 2012) and some share 
key properties of the model proposed here (Staddon et al., 1999; 
Durstewitz, 2003). An entirely different class of models is based 
on the idea that time can be read from the dynamic state of 
circuits in the cortical network (Buonomano and Mauk, 1994; 
Karmarkar and Buonomano, 2007), though the conditions for 
scalar timing in these models have not been analyzed. Some of 
the previously developed models of scalar timing are based on 
abstract entities such as counters and accumulators (Gibbon, 
1977; Church, 2003), and some are dependent on continuous 
variables (Matell and Meek, 2000; Karmarkar and Buonomano, 
2007; Oprisan and Buhusi, 2011), while others can be inter- 
preted in terms of spiking inhibitory and excitatory neurons 
(Rivest and Bengio, 2011; Simen et al., 2011) and require nearly- 
perfect integration for the decision process. The model presented 
here is formulated directly in terms populations of spiking neu- 
rons with experimentally measurable variables. There are no free 



parameters in our model, since all depend directly on neural and 
behavioral statistics. Therefore, our theory has the advantage that 
it can be tested experimentally at the physiological level. 

Our analysis indicates a very precise log-power form for the 
firing rate function. One might wonder, rightly, if it realistic to 
expect a neural processes to have such a precise formulation. It is 
important to realize that our analysis does not require or claim 
that any single neuron should display precise log-power dynam- 
ics, though to get true linear scaling the relevant population of 
neurons must possess this form. The population can be com- 
posed of individual neurons with diverse response dynamics, as 
we demonstrated in the intensity domain (Shouval et al., 2013). 
A question not answered here is how single neurons or a popu- 
lation of neurons can develop firing rate functions with a desired 
form. Possible answers are provided by previous work showing 
how single neurons with active conductances (Durstewitz, 2004; 
Shouval and Gavornik, 2011) or networks of interacting neu- 
rons (Gavornik et al., 2009; Gavornik and Shouval, 2011) can 
be tuned to, or even learn de-novo, specific temporal dynamics. 
An additional possibility is that decision neurons can select (in 
the Hebbian sense) a sub-population of existing neurons with a 
combined spike rate that has a log-power form without requir- 
ing that the dynamics of any of the individual neurons change 
at all, though we can not here propose a biologically realistic 
mechanism for making this choice. 

Recent experiments (Leon and Shadlen, 2003; Shuler and Bear, 
2006; Chubykin et al., 2013) have made physiological recordings 
from cortical cells in animals as they learn temporal discrim- 
ination tasks. These results show that the firing rate function 
of cells change when animals learn different temporal intervals 
and theoretical models have been devised to account for them 
(Reutimann et al, 2004; Gavornik et al, 2009). Analyzing these 
cases, and determining a single framework that leads to scalar 
timing, is quite different in many respects from the analysis car- 
ried out here. We are currently studying this issue. Experimentally 
it requires many trials to change the firing rate function. One 
possibility, as mentioned above, is that the model presented here 
describes mechanisms used to discriminate times in a manner 
that requires little or no learning, whereas other models would 
be required to describe how representations of specific temporal 
intervals are encoded over many trials. Regardless, the work here 
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makes a strong prediction that any neural process used to encode 
temporal intervals that display scalar timing, with our minimal 
assumptions, will have firing rates that evolve as a log-power of 
time. 
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